/*	
	This file is part of the Snoopy code.

    Snoopy code is free software: you can redistribute it and/or modify
    it under the terms of the GNU General Public License as published by
    the Free Software Foundation, either version 3 of the License, or
    (at your option) any later version.

    Snoopy code is distributed in the hope that it will be useful,
    but WITHOUT Ag.nyrtot WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
    GNU General Public License for more details.

    You should have received a copy of the GNU General Public License
    along with Snoopy code.  If not, see <http://www.gnu.org/licenses/>.
*/



#include "common.h"
#include "debug.h"

double fft_timer;

#ifdef MPI_SUPPORT

#include <fftw3-mpi.h>

fftw_plan	r2cfft_mpi_t, r2cfft_mpi, c2rfft_mpi, c2rfft_mpi_t;

void gfft_r2c_t(double *wrin) {
	double complex *win = (double complex *) wrin;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_r2c(r2cfft_mpi_t, wrin, win);
	fft_timer = fft_timer + get_c_time();
	return;
}

void gfft_c2r_t(double complex *win){
	double *wrin = (double *) win;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_c2r(c2rfft_mpi_t, win, wrin);
	fft_timer = fft_timer + get_c_time();
	return;
}

void gfft_r2c(double *wrin) {
	double complex *win = (double complex *) wrin;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_r2c(r2cfft_mpi, wrin, win);
	fft_timer = fft_timer + get_c_time();
	return;
}

void gfft_c2r(double complex *win){
	double *wrin = (double *) win;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_c2r(c2rfft_mpi, win, wrin);
	fft_timer = fft_timer + get_c_time();
	return;
}


void init_gfft() {
	
	DEBUG_START_FUNC;
	
	
#ifdef _OPENMP
	fftw_plan_with_nthreads( nthreads );
#endif
	
	r2cfft_mpi_t = fftw_mpi_plan_dft_r2c_3d( g.nyrtot, g.nxrtot, g.nzrtot, wr1, w1,  MPI_COMM_WORLD, FFT_PLANNING | FFTW_MPI_TRANSPOSED_OUT);	
	if (r2cfft_mpi_t == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW R2C_T plan creation failed");
	
	r2cfft_mpi = fftw_mpi_plan_dft_r2c_3d( g.nxrtot, g.nyrtot, g.nzrtot, wr1, w1,  MPI_COMM_WORLD, FFT_PLANNING);
	if (r2cfft_mpi == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW R2C plan creation failed");

	c2rfft_mpi_t = fftw_mpi_plan_dft_c2r_3d( g.nyrtot, g.nxrtot, g.nzrtot, w1, wr1,  MPI_COMM_WORLD, FFT_PLANNING | FFTW_MPI_TRANSPOSED_IN);
	if (c2rfft_mpi_t == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW C2R_T plan creation failed");
	
	c2rfft_mpi = fftw_mpi_plan_dft_c2r_3d( g.nxrtot, g.nyrtot, g.nzrtot, w1, wr1,  MPI_COMM_WORLD, FFT_PLANNING);
	if (c2rfft_mpi == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW C2R plan creation failed");
	 
	// init transpose routines (These are used by remap routines)
	init_transpose();
	
	fft_timer=0.0;
	
	DEBUG_END_FUNC;
	
	return;
}

void finish_gfft() {
	DEBUG_START_FUNC;
	
	fftw_destroy_plan(r2cfft_mpi_t);
	fftw_destroy_plan(r2cfft_mpi);
	fftw_destroy_plan(c2rfft_mpi_t);
	fftw_destroy_plan(c2rfft_mpi);
	
	finish_transpose();
	DEBUG_END_FUNC;
	
	return;
}


/**********************************************************************
***********************************************************************
*********     N O    M P I    R O U T I N E S        ******************
***********************************************************************
***********************************************************************/
// These routines are essentially wrappers for fftw3 routines

#else

// Here, we assume we don't have MPI

fftw_plan	r2cfft, c2rfft;

void gfft_r2c_t(double *wrin) {
	double complex *win = (double complex *) wrin;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_r2c(r2cfft, wrin, win);
	fft_timer = fft_timer + get_c_time();
	return;
}

void gfft_c2r_t(double complex *win){
	double *wrin = (double *) win;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_c2r(c2rfft, win, wrin);
	fft_timer = fft_timer + get_c_time();
	return;
}

void gfft_r2c(double *wrin) {
	double complex *win = (double complex *) wrin;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_r2c(r2cfft, wrin, win);
	fft_timer = fft_timer + get_c_time();
	return;
}

void gfft_c2r(double complex *win){
	double *wrin = (double *) win;
	fft_timer = fft_timer - get_c_time();
	fftw_execute_dft_c2r(c2rfft, win, wrin);
	fft_timer = fft_timer + get_c_time();
	return;
}

void init_gfft() {
	double complex *wi1;
	double *wir1;
	
	DEBUG_START_FUNC;
	
	wi1 = (double complex *) fftw_malloc( sizeof(double complex) * g.nloc);
	if (wi1 == NULL) ERROR_HANDLER( ERROR_CRITICAL, "No memory for wi1 allocation");
	
	wir1 = (double *) wi1;
	
#ifdef _OPENMP
	fftw_plan_with_nthreads( nthreads );
#endif
	
#ifdef WITH_2D
	r2cfft = fftw_plan_dft_r2c_2d( g.nxrtot, g.nyrtot, wr1, w1,  FFT_PLANNING);
	if (r2cfft == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW R2C plan creation failed");

	c2rfft = fftw_plan_dft_c2r_2d( g.nxrtot, g.nyrtot, w1,  wr1, FFT_PLANNING);
	if (c2rfft == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW C2R plan creation failed");
#else
	r2cfft = fftw_plan_dft_r2c_3d( g.nxrtot, g.nyrtot, g.nzrtot, wr1, w1,  FFT_PLANNING);
	if (r2cfft == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW R2C plan creation failed");

	c2rfft = fftw_plan_dft_c2r_3d( g.nxrtot, g.nyrtot, g.nzrtot, w1,  wr1, FFT_PLANNING);
	if (c2rfft == NULL) ERROR_HANDLER( ERROR_CRITICAL, "FFTW C2R plan creation failed");
#endif

	
	fftw_free(wi1);
	
	fft_timer=0.0;
	
	DEBUG_END_FUNC;
	
	return;
}

void finish_gfft() {
	DEBUG_START_FUNC;
	
	fftw_destroy_plan(r2cfft);
	fftw_destroy_plan(c2rfft);
	
	DEBUG_END_FUNC;
	
	return;
}

#endif

double read_fft_timer() {
	return(fft_timer);
}
	
